Supplementary Material
Spatial and temporal variability in the intrinsic productivity of Antarctic krill (Euphausia superba) along the Western Antarctic Peninsula under different environmental and life history scenarios
Mardones, M., Watters, G., Cárdenas, C. & Santa Cruz, F.
30 April, 2024
## Installing the Package
#The LBSPR package is now available on CRAN:
#install.packages("LBSPR")
#install.packages("devtools")
#devtools::install_github("AdrianHordyk/LBSPR")
###load the package
library(LBSPR)
library(devtools)#to install_github
library(dplyr)
library(tidyr)
library(ggplot2)
library(stringr)
library(ggpubr)
library(kableExtra)
library(hrbrthemes)
library(ggthemes)
library(tidyverse)My_theme <- theme(axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
axis.ticks.x=element_blank(),
panel.background = element_blank(),
panel.border = element_rect(fill = NA, size = 1),
strip.background = element_rect(fill = "white",
color = "white",
size = 1),
text = element_text(size = 14),
panel.grid.major = element_line(colour = "white",
size = 0.1),
panel.grid.minor = element_line(colour = "white",
size = 0.1))CONTEXT
This supplementary material with essential formulas, datasets, and code snippets is delivered for the application of the LBSPR model in studying krill populations within the West Antarctic Peninsula (WAP) region. The LBSPR model, pivotal for unraveling krill dynamics, is reinforced by mathematical expressions elucidating key parameters and environmental influences. Curated datasets provide valuable insights into factors shaping krill habitats. With a strong emphasis on reproducibility, this resource ensures that findings can undergo independent validation, fostering trust in scientific outcomes. Transparent documentation of methods and assumptions promotes understanding and scrutiny of the analysis process by peers.
METHODOLOGY
Study area
The study area includes subarea 48.1, which is one of the sectors where today the largest amount of krill fishing is concentrated (Figure 1).
Figure 1: Subarea 48.1 and management strata considered in the spatio-temporal analysis of intrinsic productivity of Krill (BS=Brainsfield Strait, EI= Elephant Island, GS= Gerlache Strait, JOIN= Joinville Island, SSWI= South West)
Monitoring Data (SISO Program)
For this analysis, data from the monitoring of the krill fishery were used, which have been systematically collected on board fishing vessels by the CCAMLR SISO (Scheme of International Scientific Observation) program. Krill sizes compositions were obtained from the entire area 48.1, which was combined in each management stratum defined at 2.1 section (Figure 2). For this analysis, 1,266,267 krill individuals were measured from fishery activity, the majority (~75%) from the Bransfield, Elephant and SouthWest strata.
Figure 2: Sizes compositions from SISO program monitoring krill fishery by strata (BS=Brainsfield Strait, EI= Elephant Island, GS= Gerlache Strait, JOIN= Joinville Island, SSWI= South West). Red line represent recruit size
The information gaps (years without sizes composition data) are not calculated because there is no autocorrelation between years, but singular estimators over time.
Biological and fishery parameters krill
#Create a new LB_pars Object
#To create a new LB_pars object you use the new function:
MyPars <- new("LB_pars")
#You can see the elements or slots of the LB_pars object using the
#slotNames function:
#slotNames(MyPars)The model needs specifications related to both biological and fishery parameters according to species evaluated. In a descriptive way, the main parameter sets are described as follows;
Biology
- von Bertalanffy asymptotic length
Linf - M/K ratio (natural mortality)divided by von Bertalanffy K coefficient)
MK - Length at 50% maturity (
L50) - Length at 95% maturity (
L95)
Fishery
- Length at 50% selectivity (
SL50) - Length at 95% selectivity (
SL95) - Biological Reference Point (BRP). F/M ratio (
FM) or Spawning Potential Ratio (SPR). If you specify both, the F/M value will be ignored.
Size Classes
-Width of the length classes (BinWidth)
MyPars <- new("LB_pars")
## A blank LB_pars object created
## Default values have been set for some parameters
MyPars@Species <- "Euphausia superba"
MyPars@Linf <- 60
MyPars@L50 <- 34
MyPars@L95 <- 55 # verrificar bibliografia
MyPars@MK <- 0.4/0.45
#Explotation
MyPars@SL50 <- 40#numeric() #1
MyPars@SL95 <- 56#numeric() #27
MyPars@SPR <- 0.75 #numeric()# ###cambia el numero 0.4 a en blanco
MyPars@BinWidth <- 1
#MyPars@FM <- 1
MyPars@Walpha <- 1
MyPars@Wbeta <- 3.0637 #r2 = 0.9651
MyPars@BinWidth <-1
MyPars@BinMax <- 70
MyPars@BinMin <- 0
MyPars@L_units <- "mm"These parameters were taken from previous works about krill life history and fishery (Maschette et al., 2020; Thanassekos et al., 2014), which are described (Table 1.
tablepar <-data.frame(Value=c(MyPars@Linf,
MyPars@L50 ,
MyPars@L95 ,
round(MyPars@MK,3) ,
MyPars@SL50 ,
MyPars@SL95 ,
MyPars@SPR ,
MyPars@Walpha ,
MyPars@Wbeta ,
MyPars@BinMax,
MyPars@BinMin ,
MyPars@L_units),
Descrption=c("VB asymptotic length",
"Maturity 50%",
"Maturity 95%",
"M/K Ratio",
"Selectivity 50%",
"Seletivity 95%",
"SPR",
"a (Length-Weight Relation)",
"b (Length-Weight Relation)",
"Bin Min",
"Bin Max",
"Units"))kbl(tablepar,
longtable = F,
booktabs = T,
caption = "\\label{Table1}Krill biological and fishery parameters") %>%
kable_styling(latex_options = c("striped",
"hold_position"),
html_font = "arial") | Value | Descrption |
|---|---|
| 60 | VB asymptotic length |
| 34 | Maturity 50% |
| 55 | Maturity 95% |
| 0.889 | M/K Ratio |
| 40 | Selectivity 50% |
| 56 | Seletivity 95% |
| 0.75 | SPR |
| 1 | a (Length-Weight Relation) |
| 3.0637 | b (Length-Weight Relation) |
| 70 | Bin Min |
| 0 | Bin Max |
| mm | Units |
.
Model Estimation LBSPR
Recent work has shown that under equilibrium conditions (that is, constant F and no recruitment variability) and assuming the von Bertalanffy growth equation, constant natural mortality for all ages, and logistic or jack-knife selectivity, standardization of the composition of lengths of two populations with the same ratio of natural mortality to growth rate (M/k) and the same ratio of mortality by fishing to natural mortality (F/M) will be identical (A. R. Hordyk et al., 2016). Extension of this model to incorporate length-at-age variability and logistic selectivity confirms that, at equilibrium, the composition of the predicted duration of catch of an exploited population is primarily determined by the ratios of M/k and F/M. The analytical models developed in A. Hordyk et al. (2014) suggest that with knowledge of the asymptotic von Bertalanffy length \(L_{\infty}\) and the coefficient of variation in \(CVL_{\infty}\), the ratio of total mortality to the von Bertalanffy growth coefficient (Z/k) for a given population can be estimated from a representative sample of the size structure of the catch. If M/k (or parameters) is also known, then the results of A. R. Hordyk et al. (2016) suggest that it is possible to estimate F/M from the composition of the catch. Often the F/M ratio has been used as a biological reference point when is 1 (Zhou et al., 2012).
The LBSPR model requires the following parameters: an estimate of the M/k ratio, \(L_{\infty}\), \(CVL_{\infty}\), and knowledge of maturity by length (maturity ogive), both set parameters know in krill. This model uses data on composition by catch sizes to estimate intrinsic productivity or Spawning Potential Ratio (SPR). This concept was extracted from Goodyear (1993), where the ratio of lifetime average egg production per recruit (EPR) was calculated for fished and no fished (virgin condition) resources. An algorithm route for the calculation of the SPR is the following;
\[{SPR}=\frac{EPR_{fished}}{EPR_{nofished}}\] where;
\[ EPR_{fished} = \sum_{a} \begin{cases} E_{a}, a = 0 \\ e^{-Z_{{a-1}}a} {E_a}, 0 < a < a_\le{max} \end{cases} \]
where \(Z_a\) = \(M+F_a\), and \(E\) is egg production at age assumed to be proportional to weight;
\[E_a \in Mat_a W_a\] on the other hand, the calculation of the reproductive potential is the same as that of those captured without F; \[ EPR_{nofished} = \sum_{a} E_ae^{-M_a} \] Assuming that egg production is proportional to the size of mature fish, relative fecundity-at-size is given by;
\[ Fec_{L,g} = Mat_{L,g} L^b \] where b is value of the exponent to reflect differente size fecunditity relationship and g is the fraction of recruits to group.
Assuming reasonable estimates of the M/K ratio, \(L_{\infty}\) (or \(CVL_{\infty}\)), size-at-maturity, the parameters F/M, \(SL_{50}\), and \(SL_{95}\) can be estimated from a representative sample of the length structure of the catch by minimizing the following multinomial negative loglikelihood function (NLL):
\[ NLL = argmin\sum_{i} O_i ln\frac{\hat{P}_i}{\hat{O}_i} \]
where \(O_i\) and \(\hat{O}_i\) are the observed number and proportion in length class i, respectively, and \(\hat{P}_i\) is the model estimate of the probability in length class i (A. R. Hordyk et al., 2016).
Like any assessment method, the LBSPR model relies on a number of simplifying assumptions. In particular, the LBSPR models are both equilibrium based and that the length composition data is representative of the exploited population at steady state (A. Hordyk et al., 2014; A. R. Hordyk et al., 2016). This methodology was implemented through the package (A. Hordyk et al., 2014) and this code and data could be visited in LBSPRKrill.
Biological References Point (PBR) in krill.
A constant challenge for krill has been to provide indicators of stock status that can be compared to predetermined biological reference points. The intrinsic productivity or SPR is commonly used to set the limit and target reference point (Goodyear, 1993; Mace, 2001; Prince et al., 2014). By definition, the SPR is equal to 100% in an unexploited stock, and zero in a non-spawning stock (eg all mature fish have been removed, or all females have been caught). The F40%, that is, the fishing mortality that allows an escape of 40% of the biomass to MSY, is the fishing mortality rate that translates into SPR = 40%, is considered a reference point to many species (Clark, 2002). Suitable SPR biological points can be derived from hypotheses about the steepness of the stock-recruit relationship (Brooks, 2013; A. R. Hordyk et al., 2016).
However Prince et al. (2014) provide some considerations in relation to the life strategy of the organisms and the SPR necessary to ensure the sustainability of the fishery. In this work, three types (I, II and III) are identified, which refers to the r and K life strategy. The krill is Type 1 (r-strategy), M/k=~1 (m=0.4, k =0.43) therefore, it requires a higher save of SPR.
On the other hand, the current krill fisheries management schemes established by CCAMLR (2010), Constable et al. (2000) proposes a decision rule scheme based on biomass that ensure the sustainability of this resource in the Southern Ocean. Through a simulation process, a 20-year projection of the krill population is generated, and a probability distribution is established around two decision rules. The first rule is based on the lowest level of exploitation allowed, around 20% of biomass, which could be consider a limit reference point. The second is the target level of the fishery, which indicates the statistical distribution of the biomass at the end of the 20-year projection under a constant catch that allows the average escape of 75% at pre-exploited levels (CCAMLR-2000 Survey). Linking this management scheme to intrinsic productivity of krill, we used two reference points for intrinsic krill productivity, 20% SPR as the limit reference point and 75% as the target to this fishery.
With this references points we can test a challenge for fishing exploitation and sustainability of krill in the Antarctic Peninsula, which is related to overfishing by recruitment. This condition occurs when the adult (reproductive) fraction has been severely reduced due to excessive fishing, jeopardizing the success of recruitment and stock replenishment in order to ensure sustainability (Cubillos, 2005; Hilborn & Walters, 1992) in relation to the previously mentioned PBRs.
The LBSPR package can be used to generate the expected size composition, the SPR, and relative yield for a given set of biological and exploitation pattern parameters. The output of the LBSPRsim function is an object of class LB_obj. This is another LBSPR object, and contains all of the information from the LB_pars object and the output of the LBSPRsim function.
It is important to note that the F/M ratio reported in the LBSPR model refers to the apical F over the adult natural mortality rate. That is, the value for fishing mortality refers to the highest level of F experienced by any single size class (A. Hordyk et al., 2014). If the selectivity pattern excludes all but the largest individuals from being exploited, it is possible to have a very high F/M ratio in a sustainable fishery (high SPR) and visceverse. Note that only the life history parameters need to be specified for the estimation model. The exploitation parameters will be estimated (A. Hordyk et al., 2014).
## Plotting the Simulation
MySim <- LBSPRsim(MyPars,
Control=list(modtype="GTG"))
# Plot with
plotSim(MySim, type="len.freq",
Cols = NULL,
size.axtex = 12,
size.title = 14,
size.SPR = 4,
size.leg = 12,
inc.pts = TRUE,
size.pt = 6)Two objects are required to fit the LBSPR model to length data: life-history parameters (LB_pars) described previously (2.4. section) and size compositions (LB_lengths), which contains the length frequency data. We provide a set of global size data and also by strata, with which we will be able to identify the spatial differences of the potential.
MyLengths <- new("LB_lengths")
slotNames(MyLengths)However, it is probably easier to create the LB_lengths object by directly reading in a .csv. A length frequency data of krill set with multiple years (2001-2020) by strata in 48.1 Subarea.
#Brainsflied Strata
datdir <- setwd("~/DOCAS/LBSPR_Krill")
Lenbs <- new("LB_lengths", LB_pars=MyPars, file=paste0(datdir, "/Length_481_Krill_2.csv"), dataType="freq",sep=";",header=T)
#Elephan Island Strata
Lenei <- new("LB_lengths", LB_pars=MyPars, file=paste0(datdir, "/lenghtEI.csv"), dataType="freq",sep=";",header=T)
#Gerlache Strata
Lengc <- new("LB_lengths", LB_pars=MyPars, file=paste0(datdir, "/lenghtGerlache.csv"), dataType="freq",sep=";",header=T)
#Join Strata
Lenjo <- new("LB_lengths", LB_pars=MyPars, file=paste0(datdir, "/lenghtJOIN.csv"), dataType="freq",sep=";",header=T)
#SSIW Strata
Lenssiw <- new("LB_lengths", LB_pars=MyPars, file=paste0(datdir, "/lenghtSSIW.csv"), dataType="freq",sep=";",header=T)## 2.9. Fit the Model by strata
# The LBSPR model is fitted by strata using the `LBSPRfit` function.
fitbs <- LBSPRfit(MyPars, Lenbs)
fitei <- LBSPRfit(MyPars, Lenei)
fitgc <- LBSPRfit(MyPars, Lengc)
fitjo <- LBSPRfit(MyPars, Lenjo)
fitssiw <- LBSPRfit(MyPars, Lenssiw)Sensitivity and perfomance analysis metohodology.
Ten sensitivity scenarios based on the upper and lower range for the asymptotic length von Bertalanffy \(L_{\infty}\) (55 to 65 mm) used in the base model (60 mm) were tested to identify the impact of this parameter on the SPR estimation, given that it is the parameters exhibiting the high degree of variability and one of the factors that most determines the estimates. On the other hand, the interdependence between krill and their environment is a well-known and influential factor in population dynamics, ecosystem impacts, and fishery. This interdependence also affects the reproductive potential and consequently, to any management decision that takes into account population parameters of krill. To assess the impact of individual growth variability, three scenarios of the VB k parameter were tested, representing different growth types (low = 0.2, medium = 0.7, and high = 1.2). Figure 3 displays a theoretical growth curve for krill based on three scenarios that were tested using LBSPR.
# Definir los parámetros del modelo de Von Bertalanffy
L_inf <- 60
k <- c(0.2, 0.7, 1.2)
t_0 <- 0
# Crear un vector de edades
edades <- seq(0, 8, by = 0.2)
# Calcular los valores teóricos de crecimiento somático
df <- data.frame(edades = rep(edades,
length(k)),
k = rep(k, each = length(edades)))
df <- df %>%
mutate(crecimiento = L_inf * (1 - exp(-k * (edades - t_0))))
# Graficar la curva de crecimiento somático
sen <- ggplot(df, aes(x = edades,
y = crecimiento,
color = factor(k))) +
geom_line(size=1.5) +
xlab("Ages") +
ylab("Length (mm)") +
ggtitle("") +
theme(panel.background = element_rect(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())+
scale_color_viridis_d(option="C",
name = "Theoretical growth",
labels = c("Low", "Medium", "High"))+
theme_few()
sen
Figure 3: Theoretical growth scenarios to krill varability
After applying each scenario using each of the parameter setting, the results of scenarios are compared with the results provided by the methods based setting (Table 1), analyzing in this way the efect of underestimation/overestimation of the parameters \(L_{\infty}\) and k.
RESULTS
Model Perfomance
The results demonstrate good fits of the size structure distribution model for krill across the strata. The model accurately captures the distribution patterns of size classes, indicating its effectiveness in characterizing the population structure. The model successfully captures the variations in size distribution, reflecting the natural variability in krill populations across different strata. These findings suggest that the model can be utilized as a valuable tool for understanding and predicting size structure dynamics in krill populations (Figure 4).
plotSize(fitei,
Title="Braisnfield strata")+
My_theme
Figure 4: Fit of the model to the data of lengths in Braishfield strata
The Joinville strata exhibits the least model fits, supposedly attributed to the limited available size compositions data spanning only three years (Figure 5).
plotSize(fitjo,
Title="Joinville strata")+
My_theme
Figure 5: Fit of the model to the data of lengths in Joinville strata
For all the other strata (Gerlache, EI and SSWI) the model predict sizes compositions in a correct way for all the years (for details see LBSPRKrill).
plotSize(fitei,
Title="Elephan Island strata")
plotSize(fitgc,
Title="Gerlache strata")
plotSize(fitjo,
Title="Joinville strata")
plotSize(fitssiw,
Title="SSWI strata")The difference between the observed accumulated size compositions for each stratum and compare it with the expected size composition at a target SPR (75% SPR). In the simulation of the structure in its virgin condition (without fishing), the red bars represent each stratum. Additionally, the overlap with the average structures observed during the years of fishery monitoring can be visualized. The SSWI Gerlache and EI strata exhibit the greatest differences from the simulated structure, possibly due to the significant contribution of juveniles in these strata. Conversely, the BS and JO strata demonstrate the closest resemblance to the simulated structure (Figure 6).
yr <- 5 # first year of data
bscom <- plotTarg(MyPars, Lenbs, yr=yr,
Cols = c(2,4),
size.axtex = 8,
title="BS",
targtext = FALSE)+
My_theme
eicom <- plotTarg(MyPars, Lenei, yr=yr,
Cols = c(2,4),
size.axtex = 8,
title="EI",
targtext = FALSE)+
My_theme
gccom <- plotTarg(MyPars, Lengc, yr=yr,
Cols = c(2,4),
size.axtex = 8,
title="GS",
targtext = FALSE)+
My_theme
jocom <- plotTarg(MyPars, Lenjo, yr=1,
Cols = c(2,4),
size.axtex = 8,
title="JOIN",
targtext = FALSE)+
My_theme
sscom <- plotTarg(MyPars, Lenssiw , yr=yr,
Cols = c(2,4),
size.axtex = 8,
title="SSWI",
targtext = FALSE)+
My_theme
ggarrange(bscom, eicom + rremove("ylab"),
gccom + rremove("ylab"),
jocom, sscom + rremove("ylab"),
ncol = 3, nrow = 2,
legend="bottom",
common.legend = TRUE)
Figure 6: Difference between the observed accumulated size structure for each stratum related SPR objective
Furthermore, ogive maturity curve specific to each stratum and year, as well as the estimated length selectivity curve, are presented in (Figure 7). These curves offer crucial insights into the fishery’s impact on the population and its reproductive status, providing measures of the population’s vulnerability to fishing mortality. Notably, the Braisnfield stratum stands out with a lower proportion of mature individuals, suggesting a higher prevalence of juveniles. It is important to note that the same maturity parameters were applied across all strata.
bsmat <- plotMat(fitbs,
useSmooth = TRUE,
Title="BS")+
My_theme
eimat <- plotMat(fitei,
useSmooth = TRUE,
Title="EI")+
My_theme
gcmat <- plotMat(fitgc,
useSmooth = TRUE,
Title="GS")+
My_theme
jomat <- plotMat(fitjo,
useSmooth = TRUE,
Title="JOIN")+
My_theme
ssmat <- plotMat(fitssiw,
useSmooth = TRUE,
Title="SSWI")+
My_theme
ggarrange(bsmat,
eimat + rremove("ylab"),
gcmat ,
jomat + rremove("ylab"),
ssmat,
ncol = 2, nrow = 3,
legend="right",
common.legend = TRUE)
Figure 7: Maturity curves by strata
Comparing producivity between years and Stratas
The analysis of the krill population’s reproductive potential across different years and strata reveals significant differences. Brainsfield and Gerlache strata exhibit a low reproductive potential below the proposed management target of 75% in the last year, with values of 0.121 and 0.085, respectively, falling even below the limit reference point. This condition arises from the concentration of a substantial number of immature individuals (juveniles) in these strata, which are being exploited by the fishery, thereby hindering their reproductive cycles from completing. On the contrary, the Elephant Island stratum demonstrates higher spawning potential ratio (SPR) levels in recent years, reaching 0.421 in 2019, which aligns closer to the management objective. This discrepancy is attributed to the spatial distribution of krill, as the Elephant Island stratum possesses a larger proportion of adult individuals compared to other strata. Figure 8 provides a visual representation of the SPR trends across years and strata, clearly indicating the references (yellow line = 75% SPR Objective and Red line = 20% Limit SPR).
# Get values estimated
sprbs <- as.data.frame(cbind(fitbs@Years,
fitbs@SPR,
fitbs@Vars[,4]))
colnames(sprbs) <- c("Year","SPR", "Var")
sprbs$SPRv <- rep("BS", nrow(sprbs))
sprei <- as.data.frame(cbind(fitei@Years,
fitei@SPR,
fitei@Vars[,4]))
colnames(sprei) <- c("Year","SPR", "Var")
sprei$SPRv <- rep("EI", nrow(sprei))
sprgc <- as.data.frame(cbind(fitgc@Years,
fitgc@SPR,
fitgc@Vars[,4]))
colnames(sprgc) <- c("Year","SPR", "Var")
sprgc$SPRv <- rep("GS", nrow(sprgc))
sprjo <- as.data.frame(cbind(fitjo@Years,
fitjo@SPR,
fitjo@Vars[,4]))
colnames(sprjo) <- c("Year","SPR", "Var")
sprjo$SPRv <- rep("JOIN", nrow(sprjo))
sprssiw <- as.data.frame(cbind(fitssiw@Years,
fitssiw@SPR,
fitssiw@Vars[,4]))
colnames(sprssiw) <- c("Year","SPR","Var")
sprssiw$SPRv <- rep("SSWI", nrow(sprssiw))
allspr <- rbind(sprbs, sprei, sprgc, sprjo, sprssiw)
allsprwide <- round(pivot_wider(allspr,
names_from = SPRv,
values_from = c(SPR, Var),
values_fill = NA,
names_sep = " ",),3)allsprpl <- ggplot(allspr,
aes(Year,
SPR,
color=SPRv))+
geom_point(alpha=3)+
stat_smooth(method = "lm",
alpha=0.3,
se=FALSE)+
geom_hline(yintercept = 0.75,
colour= '#006d2c',
alpha=0.5,
linewidth=1.2)+
geom_hline(yintercept = 0.20,
colour= '#bd0026',
alpha=0.5,
linewidth=1.2)+
facet_wrap(.~SPRv, ncol = 5)+
scale_color_viridis_d(option = "F")+
xlim(2001,2021)+
ylim(0,1.2)+
theme_minimal()+
theme(legend.position = "none",
axis.text.x = element_text(angle = 90, hjust = 2),
panel.background = element_rect(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())
allsprpl
Figure 8: Krill Intrinsic Productivity (SPR) by strata and by year
The Elephant Island stratum has exhibited a higher prevalence of adult fraction in fishing catches, leading to an increase in reproductive potential in recent years. Conversely, the Bransfield and Gerlache strata experience intense recruitment overfishing, with their reproductive potential falling significantly below the recommended target at 54% and 66.5%, respectively.
All the estimated values of SPR and their associated variance by stratum and by year can be identified in Table 2.
kbl(allsprwide,
longtable = F,
booktabs = T,
caption = "\\label{Table2}Estimates SPR by Strata") %>%
add_header_above(c(" ", "SPR" = 5, "Variance" = 5)) %>%
kable_styling(latex_options = c("scale_down", "hold_position"))| Year | SPR BS | SPR EI | SPR GS | SPR JOIN | SPR SSWI | Var BS | Var EI | Var GS | Var JOIN | Var SSWI |
|---|---|---|---|---|---|---|---|---|---|---|
| 2001 | 0.339 | 0.223 | NA | NA | 0.194 | 0.025 | 0.017 | NA | NA | 0.003 |
| 2004 | 0.209 | 0.221 | NA | NA | 0.210 | 0.003 | 0.001 | NA | NA | 0.002 |
| 2005 | 0.325 | NA | NA | NA | 0.254 | 0.001 | NA | NA | NA | 0.000 |
| 2006 | 0.314 | NA | 0.152 | NA | 0.353 | 0.026 | NA | 0.008 | NA | 0.033 |
| 2007 | 0.372 | 0.090 | 0.069 | NA | 0.274 | 0.024 | 0.002 | 0.000 | NA | 0.017 |
| 2008 | 0.125 | 0.072 | NA | NA | NA | 0.002 | 0.001 | NA | NA | NA |
| 2009 | 0.296 | 0.319 | NA | NA | 0.234 | 0.003 | 0.011 | NA | NA | 0.001 |
| 2010 | 0.310 | 0.424 | 0.246 | NA | 0.457 | 0.006 | 0.027 | 0.019 | NA | 0.013 |
| 2011 | 0.579 | NA | NA | NA | 0.417 | 0.010 | NA | NA | NA | 0.014 |
| 2012 | 0.395 | 0.176 | 0.149 | NA | 0.362 | 0.061 | 0.006 | 0.001 | NA | 0.017 |
| 2013 | 0.193 | 0.182 | 0.128 | NA | 0.182 | 0.001 | 0.001 | 0.000 | NA | 0.003 |
| 2014 | 0.235 | 0.157 | 0.177 | NA | 0.294 | 0.001 | 0.000 | 0.001 | NA | 0.002 |
| 2015 | 0.247 | 0.385 | 0.182 | NA | 0.378 | 0.005 | 0.016 | 0.001 | NA | 0.023 |
| 2016 | 0.331 | 0.317 | 0.294 | 0.331 | 0.301 | 0.029 | 0.029 | 0.003 | 0.036 | 0.026 |
| 2017 | 0.211 | 1.000 | 0.190 | 0.254 | 0.262 | 0.008 | 0.000 | 0.015 | 0.002 | 0.007 |
| 2018 | 0.216 | 0.330 | 0.132 | NA | 0.291 | 0.002 | 0.034 | 0.001 | NA | 0.033 |
| 2019 | 0.331 | 0.421 | 0.125 | 0.345 | 0.367 | 0.010 | 0.018 | 0.002 | 0.023 | 0.002 |
| 2020 | 0.121 | NA | 0.085 | NA | 0.240 | 0.002 | NA | 0.001 | NA | 0.001 |
Sensitivity and perfomance analysis
The results of the methods in the reference setting are compared to the obtained under overstimation/underestimation in intrinsic productivity regarding asymptotic length von Bertalanffy \(L_{\infty}\) parameter in krill. First, it was possible to identify that for all strata, the impact of low \(L_{\infty}\) ranges under the base model (60 mm) overestimates the level of reproductive potential krill with values between 42% and 32% (Bransfield and Elephant Island strata respectively). Regarding higher \(L_{\infty}\) settings, the model tends to underestimate the reproductive potential with values between -25% and -30% (Gerlache and Joinville strata) Figure 9, Table 3.
for (i in 55:65) {
nombre_objeto <- paste0("MyPars_", i) # Crear un nombre único para cada objeto modificado
objeto_modificado <- MyPars # Crear una copia del objeto original
objeto_modificado@Linf <- i # Modificar el valor de @Linf en cada iteración
assign(nombre_objeto, objeto_modificado, envir = .GlobalEnv) # Asignar el objeto modificado al nombre creado
# Resto del código que deseas ejecutar con el objeto modificado
## A blank LB_pars object created
## Default values have been set for some parameters
MyPars@Species <- "Euphausia superba"
MyPars@L50 <- 34
MyPars@L95 <- 55 # verrificar bibliografia
MyPars@MK <- 0.4/0.45
#Explotation
MyPars@SL50 <- 40#numeric() #1
MyPars@SL95 <- 56#numeric() #27
MyPars@SPR <- 0.75 #numeric()# ###cambia el numero 0.4 a en blanco
MyPars@BinWidth <- 1
#MyPars@FM <- 1
MyPars@Walpha <- 1
MyPars@Wbeta <- 3.0637 #r2 = 0.9651
MyPars@BinWidth <-1
MyPars@BinMax <- 70
MyPars@BinMin <- 0
MyPars@L_units <- "mm"
# Imprimir el valor de @Linf después de cada cambio
#print(get(nombre_objeto)@Linf)
}# Recuerdo las bases de datos de lengths
#Brainsflied Strata: Lenbs
#Elephan Island Strata: Lenei
#Gerlache Strata:Lenex
#Join Strata: Lenjo
#SSIW Strata: Lenssiw
# ciclo de ajuste para BS
for (i in 55:65) {
nombre_objeto <- paste0("MyPars_", i) # Nombre del objeto S4 en cada iteración
nombre_variable <- paste0("fitbs", i) # Nuevo nombre para el resultado
# Obtener el objeto S4 correspondiente
objeto <- get(nombre_objeto)
# Ejecutar LBSPRfit() y asignar el resultado a una nueva variable con nombre único
assign(nombre_variable, LBSPRfit(objeto, Lenbs), envir = .GlobalEnv)
}
# ciclo de ajuste para Gerlache
for (i in 55:65) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitei", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lenei), envir = .GlobalEnv)
}
# ciclo de ajuste para Gerlache
for (i in 55:65) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitgc", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lengc), envir = .GlobalEnv)
}
# ciclo de ajuste para JoinVIlle
for (i in 55:65) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitjo", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lenjo), envir = .GlobalEnv)
}
# ciclo de ajuste para JSSWI
for (i in 55:65) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitsswi", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lenssiw), envir = .GlobalEnv)
}#Genero lo que quiero comparar BS
valbs <- as.data.frame(cbind(fitbs55@Years,
fitbs55@SPR,
fitbs56@SPR,
fitbs57@SPR,
fitbs58@SPR,
fitbs59@SPR,
fitbs60@SPR,
fitbs61@SPR,
fitbs62@SPR,
fitbs63@SPR,
fitbs64@SPR,
fitbs65@SPR))
colnames(valbs) <- c("Years", "Linf55","Linf56", "Linf57", "Linf58",
"Linf59", "Linf60",
"Linf61", "Linf62" , "Linf63", "Linf64",
"Linf65")
valbs_largo <- pivot_longer(valbs,
cols = c(2:12,),
names_to = "Parameter",
values_to = "SPR")
valbs_largo$SPRstra <- rep("BS", nrow(valbs_largo))
#Genero lo que quiero comparar EI
valei <- as.data.frame(cbind(fitei55@Years,
fitei55@SPR,
fitei56@SPR,
fitei57@SPR,
fitei58@SPR,
fitei59@SPR,
fitei60@SPR,
fitei61@SPR,
fitei62@SPR,
fitei63@SPR,
fitei64@SPR,
fitei65@SPR))
colnames(valei) <- c("Years", "Linf55","Linf56", "Linf57", "Linf58",
"Linf59", "Linf60",
"Linf61", "Linf62" , "Linf63", "Linf64",
"Linf65")
valei_largo <- pivot_longer(valei,
cols = c(2:12,),
names_to = "Parameter",
values_to = "SPR")
valei_largo$SPRstra <- rep("EI", nrow(valei_largo))
#Genero lo que quiero comparar Gerlache
valgc <- as.data.frame(cbind(fitgc55@Years,
fitgc55@SPR,
fitgc56@SPR,
fitgc57@SPR,
fitgc58@SPR,
fitgc59@SPR,
fitgc60@SPR,
fitgc61@SPR,
fitgc62@SPR,
fitgc63@SPR,
fitgc64@SPR,
fitgc65@SPR))
colnames(valgc) <- c("Years", "Linf55","Linf56", "Linf57", "Linf58",
"Linf59", "Linf60",
"Linf61", "Linf62" , "Linf63", "Linf64",
"Linf65")
valgc_largo <- pivot_longer(valgc,
cols = c(2:12,),
names_to = "Parameter",
values_to = "SPR")
valgc_largo$SPRstra <- rep("GS", nrow(valgc_largo))
#Genero lo que quiero comparar JOin
valjo <- as.data.frame(cbind(fitjo55@Years,
fitjo55@SPR,
fitjo56@SPR,
fitjo57@SPR,
fitjo58@SPR,
fitjo59@SPR,
fitjo60@SPR,
fitjo61@SPR,
fitjo62@SPR,
fitjo63@SPR,
fitjo64@SPR,
fitjo65@SPR))
colnames(valjo) <- c("Years", "Linf55","Linf56", "Linf57", "Linf58",
"Linf59", "Linf60",
"Linf61", "Linf62" , "Linf63", "Linf64",
"Linf65")
valjo_largo <- pivot_longer(valjo,
cols = c(2:12,),
names_to = "Parameter",
values_to = "SPR")
valjo_largo$SPRstra <- rep("JOIN", nrow(valjo_largo))
#Genero lo que quiero comparar SSWI
valsswi <- as.data.frame(cbind(fitsswi55@Years,
fitsswi55@SPR,
fitsswi56@SPR,
fitsswi57@SPR,
fitsswi58@SPR,
fitsswi59@SPR,
fitsswi60@SPR,
fitsswi61@SPR,
fitsswi62@SPR,
fitsswi63@SPR,
fitsswi64@SPR,
fitsswi65@SPR))
colnames(valsswi) <- c("Years", "Linf55","Linf56", "Linf57", "Linf58",
"Linf59", "Linf60",
"Linf61", "Linf62" , "Linf63", "Linf64",
"Linf65")
valsswi_largo <- pivot_longer(valsswi,
cols = c(2:12,),
names_to = "Parameter",
values_to = "SPR")
valsswi_largo$SPRstra <- rep("SSWI", nrow(valsswi_largo))
#agrupo
valtodo <- rbind(valsswi_largo,
valjo_largo,
valbs_largo,
valei_largo,
valgc_largo)# plot all Linf
sensto <- ggplot(valtodo %>%
drop_na(),
aes(x = Parameter,
y = SPR,
fill=SPRstra)) +
geom_violin(trim=FALSE)+
geom_jitter(height = 0,
width = 0.3,
alpha=0.4,
color="grey")+
scale_fill_viridis_d(option="F")+
facet_wrap(~SPRstra,
ncol=5)+
geom_vline(xintercept = 6,color = "black")+
geom_hline(yintercept = 0.75,
colour= '#006d2c',
alpha=0.5,
linetype=2,
linewidth=1.2)+
geom_hline(yintercept = 0.20,
colour= '#bd0026',
alpha=0.5,
linetype=2,
linewidth=1.2)+
stat_summary(fun.y=mean, geom="point", size=1.1, color="red")+
labs(y="SPR Mean",
x="VB Parameters tested")+
theme_minimal()+
theme(legend.position = "none",
axis.text.x = element_text(angle = 90, hjust = 2),
panel.background = element_rect(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())+
coord_flip()
sensto
Figure 9: Sensitivity analysis by strata about asymptotic length VB
# Tables
tablinf <- valtodo %>%
group_by(Parameter, SPRstra) %>%
summarise(Median = mean(SPR),
Variance=var(SPR)) %>%
pivot_wider(names_from= SPRstra,
values_from = c(Median, Variance),
names_sep = " ") %>%
rename( "VB scenario"="Parameter") %>%
mutate_if(is.numeric, round, 2)kbl(tablinf,
longtable = F,
booktabs = T,
caption = "Estimated by asymptotyc lenght (VB) scenario") %>%
add_header_above(c(" ", "Mean" = 5, "Variance" = 5)) %>%
kable_styling(latex_options = c("scale_down", "hold_position"),
font_size = 10)| VB scenario | Median BS | Median EI | Median GS | Median JOIN | Median SSWI | Variance BS | Variance EI | Variance GS | Variance JOIN | Variance SSWI |
|---|---|---|---|---|---|---|---|---|---|---|
| Linf55 | 0.50 | 0.45 | 0.26 | 0.50 | 0.47 | 0.05 | 0.06 | 0.01 | 0 | 0.02 |
| Linf56 | 0.43 | 0.41 | 0.23 | 0.45 | 0.43 | 0.04 | 0.06 | 0.01 | 0 | 0.01 |
| Linf57 | 0.39 | 0.38 | 0.21 | 0.40 | 0.39 | 0.03 | 0.05 | 0.01 | 0 | 0.01 |
| Linf58 | 0.35 | 0.35 | 0.19 | 0.37 | 0.35 | 0.02 | 0.05 | 0.01 | 0 | 0.01 |
| Linf59 | 0.31 | 0.33 | 0.18 | 0.34 | 0.32 | 0.01 | 0.05 | 0.00 | 0 | 0.01 |
| Linf60 | 0.29 | 0.31 | 0.16 | 0.31 | 0.30 | 0.01 | 0.05 | 0.00 | 0 | 0.01 |
| Linf61 | 0.26 | 0.29 | 0.15 | 0.29 | 0.28 | 0.01 | 0.05 | 0.00 | 0 | 0.01 |
| Linf62 | 0.25 | 0.27 | 0.14 | 0.27 | 0.26 | 0.01 | 0.04 | 0.00 | 0 | 0.00 |
| Linf63 | 0.23 | 0.25 | 0.13 | 0.25 | 0.24 | 0.01 | 0.03 | 0.00 | 0 | 0.00 |
| Linf64 | 0.22 | 0.23 | 0.12 | 0.24 | 0.23 | 0.01 | 0.03 | 0.00 | 0 | 0.00 |
| Linf65 | 0.21 | 0.22 | 0.12 | 0.22 | 0.22 | 0.01 | 0.02 | 0.00 | 0 | 0.00 |
Regarding the three growth levels of the species (high, medium, low) referred to the parameter \(k\), it was possible to identify that high and medium growth types result in very low SPR (spawning potential ratio) estimates compared to slow and medium growth. In fact, with high individual growth, the model estimates that SPR levels would be very close to the target management levels (75% SPR) and far from the limit reference level of 20% (Figure 10, Table 4).
for (i in c(2,0.5,0.3)) {
nombre_objeto <- paste0("MyPars_", i) # Crear un nombre único para cada objeto modificado
objeto_modificado <- MyPars # Crear una copia del objeto original
objeto_modificado@MK <- i # Modificar el valor de @Linf en cada iteración
assign(nombre_objeto, objeto_modificado, envir = .GlobalEnv) # Asignar el objeto modificado al nombre creado
# Resto del código que deseas ejecutar con el objeto modificado
## A blank LB_pars object created
## Default values have been set for some parameters
MyPars@Species <- "Euphausia superba"
MyPars@L50 <- 34
MyPars@L95 <- 55 # verrificar bibliografia
MyPars@Linf <- 60
#Explotation
MyPars@SL50 <- 40#numeric() #1
MyPars@SL95 <- 56#numeric() #27
MyPars@SPR <- 0.75 #numeric()# ###cambia el numero 0.4 a en blanco
MyPars@BinWidth <- 1
#MyPars@FM <- 1
MyPars@Walpha <- 1
MyPars@Wbeta <- 3.0637 #r2 = 0.9651
MyPars@BinWidth <-1
MyPars@BinMax <- 70
MyPars@BinMin <- 0
MyPars@L_units <- "mm"
# Imprimir el valor de @Linf después de cada cambio
#print(get(nombre_objeto)@Linf)
}# Recuerdo las bases de datos de lengths
#Brainsflied Strata: Lenbs
#Elephan Island Strata: Lenei
#Gerlache Strata:Lengc
#Join Strata: Lenjo
#SSIW Strata: Lenssiw
# ciclo de ajuste para BS
for (i in c(2,0.5,0.3)) {
nombre_objeto <- paste0("MyPars_", i) # Nombre del objeto S4 en cada iteración
nombre_variable <- paste0("fitbs", i) # Nuevo nombre para el resultado
# Obtener el objeto S4 correspondiente
objeto <- get(nombre_objeto)
# Ejecutar LBSPRfit() y asignar el resultado a una nueva variable con nombre único
assign(nombre_variable, LBSPRfit(objeto, Lenbs), envir = .GlobalEnv)
}
# ciclo de ajuste para Gerlache
for (i in c(2,0.5,0.3)) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitei", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lenei), envir = .GlobalEnv)
}
# ciclo de ajuste para Gerlache
for (i in c(2,0.5,0.3)) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitgc", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lengc), envir = .GlobalEnv)
}
# ciclo de ajuste para JoinVIlle
for (i in c(2,0.5,0.3)) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitjo", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lenjo), envir = .GlobalEnv)
}
# ciclo de ajuste para JSSWI
for (i in c(2,0.5,0.3)) {
nombre_objeto <- paste0("MyPars_", i)
nombre_variable <- paste0("fitsswi", i)
objeto <- get(nombre_objeto)
assign(nombre_variable, LBSPRfit(objeto, Lenssiw), envir = .GlobalEnv)
}#Genero lo que quiero comparar BS
valprobs <- as.data.frame(cbind(fitbs0.3@Years,
fitbs0.3@SPR,
fitbs0.5@SPR,
fitbs2@SPR))
colnames(valprobs) <- c("Years", "High","Med", "Low")
valprobs_largo <- pivot_longer(valprobs,
cols = c(2:4,),
names_to = "Parameter",
values_to = "SPR")
valprobs_largo$SPRstra <- rep("BS", nrow(valprobs_largo))
#Genero lo que quiero comparar EI
valproei <- as.data.frame(cbind(fitei0.3@Years,
fitei0.3@SPR,
fitei0.5@SPR,
fitei2@SPR))
colnames(valproei) <- c("Years", "High","Med", "Low")
valproei_largo <- pivot_longer(valproei,
cols = c(2:4,),
names_to = "Parameter",
values_to = "SPR")
valproei_largo$SPRstra <- rep("EI", nrow(valproei_largo))
#Genero lo que quiero comparar Gerlache
valprogc <- as.data.frame(cbind(fitgc0.3@Years,
fitgc0.3@SPR,
fitgc0.5@SPR,
fitgc2@SPR))
colnames(valprogc) <- c("Years", "High","Med", "Low")
valprogc_largo <- pivot_longer(valprogc,
cols = c(2:4,),
names_to = "Parameter",
values_to = "SPR")
valprogc_largo$SPRstra <- rep("GS", nrow(valprogc_largo))
#Genero lo que quiero comparar JOin
valprojo <- as.data.frame(cbind(fitjo0.3@Years,
fitjo0.3@SPR,
fitjo0.5@SPR,
fitjo2@SPR))
colnames(valprojo) <- c("Years", "High","Med", "Low")
valprojo_largo <- pivot_longer(valprojo,
cols = c(2:4,),
names_to = "Parameter",
values_to = "SPR")
valprojo_largo$SPRstra <- rep("JOIN", nrow(valprojo_largo))
#Genero lo que quiero comparar SSWI
valprosswi <- as.data.frame(cbind(fitsswi0.3@Years,
fitsswi0.3@SPR,
fitsswi0.5@SPR,
fitsswi2@SPR))
colnames(valprosswi) <- c("Years", "High","Med", "Low")
valprosswi_largo <- pivot_longer(valprosswi,
cols = c(2:4,),
names_to = "Parameter",
values_to = "SPR")
valprosswi_largo$SPRstra <- rep("SSWI", nrow(valprosswi_largo))
#agrupo
valprotodo <- rbind(valprosswi_largo,
valproei_largo,
valprobs_largo,
valprojo_largo,
valprogc_largo)#plot
# todo
valprotodo$Parameter <- factor(valprotodo$Parameter,levels = c("Low",
"Med",
"High"))
sensproto <- ggplot(valprotodo %>%
drop_na(),
aes(x = Parameter,
y = SPR,
fill=SPRstra)) +
stat_ydensity(trim=TRUE,
position=position_dodge(1))+
geom_hline(yintercept = 0.75,
colour= '#006d2c',
alpha=0.5,
linetype="dashed")+
geom_hline(yintercept = 0.20,
colour= '#bd0026',
alpha=0.5,
linetype="dashed")+
#geom_jitter(alpha=0.4,
# color="black")+
#facet_wrap(~SPRstra,
# ncol=5)+
labs(y="SPR Mean",
x="Growth tested")+
theme_minimal()+
theme(legend.position="top",
panel.background = element_rect(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())+
scale_fill_viridis_d(option="F",
name="Strata")
sensproto
Figure 10: Sensitivity analysis by strata about krill growth type
tablk <- valprotodo %>%
group_by(Parameter, SPRstra) %>%
summarise(Median = mean(SPR),
Variance=var(SPR)) %>%
pivot_wider(names_from= SPRstra,
values_from = c(Median, Variance),
names_sep = " ") %>%
rename( "Growth scenario"="Parameter") %>%
mutate_if(is.numeric, round, 2)
kbl(tablk,
longtable = F,
booktabs = T,
caption = "Estimated by growth type scenario") %>%
add_header_above(c(" ", "Mean" = 5, "Variance" = 5)) %>%
kable_styling(latex_options = c("scale_down", "hold_position"),
font_size = 10)| Growth scenario | Median BS | Median EI | Median GS | Median JOIN | Median SSWI | Variance BS | Variance EI | Variance GS | Variance JOIN | Variance SSWI |
|---|---|---|---|---|---|---|---|---|---|---|
| Low | 0.64 | 0.60 | 0.47 | 0.71 | 0.64 | 0.02 | 0.04 | 0.03 | 0 | 0.01 |
| Med | 0.15 | 0.17 | 0.08 | 0.16 | 0.16 | 0.01 | 0.02 | 0.00 | 0 | 0.00 |
| High | 0.09 | 0.10 | 0.04 | 0.09 | 0.09 | 0.00 | 0.01 | 0.00 | 0 | 0.00 |
CODE REPOSITORY
The data, codes and another documents of this exercise can be found in the following link LBSPR-Krill